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Abstract: 

We present an exact enumeration algorithm for identifying the native configuration - a maximally 
compact self avoiding walk configuration that is also the minimum energy configuration for a given 
set of contact-energy schemes; the process is implicitly sequence-dependent. In particular, we show 
that the 25-step native configuration on a diamond lattice consists of two sheet-like structures and is 
the same for all the contact-energy schemes, {( — 1, 0, 0); (—7, —3, 0); (—7, —3, —1); (—7, —3, 1)}; 
on a square lattice also, the 24-step native configuration is independent of the energy schemes consid- 
ered. However, the designing sequence for the diamond lattice walk depends on the energy schemes 
used whereas that for the square lattice walk does not. We have calculated the temperature-dependent 
specific heat for these designed sequences and the four energy schemes using the exact density of 
states. These data show that the energy scheme (—7, —3, —1) is preferable to the other three for both 
diamond and square lattice because the associated sequences give rise to a sharp low-temperature 
peak. We have also presented data for shorter (23-, 21- and 17-step) walks on a diamond lattice to 
show that this algorithm helps identify a unique minimum energy configuration by suitably taking 
care of the ground-state degeneracy. Interestingly, all these shorter target configurations also show 
sheet-like secondary structures. 

Keywords: Lattice proteins. Self Avoiding Walks, HP-design, Specific heat 
Introduction 

Brief review - Proteins are complex biopolymers that are made up of various sequences of amino 
acids; their biological functions are intimately related to the unique and stable conformational struc- 
tures, referred to as their 'native' states. They are stable not only with respect to moderate fluc- 
tuations in environmental conditions but also with respect to minor mutations of their amino acid 
sequences jTt]. Their stability and uniqueness imply that they should correspond to deep funnel-like 
minima of the free energy landscape |12l]. Since the energy of a conformation depends on the amino 
acid sequence, it is reasonable to expect that nature favors only those sequences whose associated 
energy landscape has a pronounced global minimum. This suggests an imphcit relationship between 
the conformational state of a protein and the amino acid sequence that associates energy with that 
structure. 

There are two aspects to a study of the relationship between structure and sequence - (a) 'Folding', 
the process by which a given chain of amino acids assumes a specific native conformation of a protein 
Istl and (b) 'Design', the choice of an amino acid sequence for which a given conformation will be 
the native state of a protein . 

The number of native structures that correspond to the functional proteins is known to be much 
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less, at least by two orders of magnitude, than the total number of proteins 151]. This implies a many- 
to-one mapping of amino acid sequences onto the set of all native structures. Given the fact that an 
arbitrarily chosen sequence of amino acids cannot represent a protein, identifying the unique set of 
representative sequences that correspond to a given native structure is a challenging problem for the 
following reason. 

Since there are 20 different types of amino acids, the total number of possible sequences for a 
chain consisting of N amino acid residues is 20^; also, there are a large number of conformations, 
other than the native structure, available for a given chain of amino acid residues. Searching for a 
funnel-like global free energy minimum in the combined sequence-structure space is a formidable 
task. This necessitates building simplified models of protein that can throw light on some basic 
aspects of sequence-structure relationship. 

The observation that the hydrophobic residues are normally sequestered to the interior of a pro- 
tein's native structure implies the dominant role of hydrophobic effects in the folding process; this, 
in turn, suggests a grouping of the twenty naturally occurring amino acids into two broad classes - 
hydrophobic (H) and polar (P) This enormously reduces the size of the sequence space to be 
explored. 

By further restricting the dihedral angle to assume only a certain number of discrete values in 
the range to 27r, we ensure that the number of possible conformations of a finite chain of amino 
acids is denumerably finite. An off-shoot of such coarse-graining schemes is a lattice-protein model 
101 in which Self Avoiding Walks (SAW 1^) on a lattice are taken to represent the various possible 
conformations of an amino acid chain while a two lettered HP-sequence represents the sequence of 
residues. 

Within the framework of such a lattice-protein model, it is possible to search for a unique set of 
HP-sequences that correspond to a given SAW configuration, referred to as the 'target' configuration. 
The choice of the target configuration depends on the actual systems to be modeled and studied; for 
example. If sequence-design for compact, globular proteins is the problem of interest, then a lattice 
SAW having a maximum number of non-bonded nearest neighbor pairs, called 'contacts' could be 
chosen as the target configuration [9J. The designability of a target configuration depends on the 
nature of interactions between monomers in contact ifioll . 

Given a scheme for associating energy with a contact, it will be of general interest to devise robust 
algorithmic rules for identifying a unique target configuration along with its designing sequence. Lot 
of work reported in the literature focus on designing HP-sequence for SAW configurations within 
regular lattice shapes such as the 6x6 square or the 3x3x3 cube, and also within different shapes 
on a square lattice [1 

Different secondary structural elements are known to be favored by different lattices lfl2l ; for 
example, at low resolution (3-4 A rms), a cubic lattice favors a-helix structure whereas a diamond 
lattice favors /3-sheet structure. On the basis of a detailed comparative study of various lattices, 
Godzik et al. 11211 have shown that SAWs on lattices with higher coordination number provide 
better structural representation of a protein, even though lattice effects are minimal as far as folding 
pathways are concerned. 

Motivation - Since the various contact-energy schemes considered in the literature are based on 
statistical estimates of the average effective inter-residue interactions of real proteins lIlSll . it is of 
interest to study how sensitive is the choice of a target configuration to changes in the contact- 
energy scheme. As this depends also on the specific sequence used, choosing a target configuration 
itself should be a sequence-dependent process. In fact, the many-to-one mapping of amino acid 
sequences onto the set of all native structures implies that a given native structure should not be very 
sensitive to the set of contact-energy schemes underlying the set of sequences designing it. There 
is also a related question of realizing a funnel-shaped energy- landscape which, in a milder form, 
implies reducing the diversity of target configurations as much as possible. Subsequently, for a target 
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configuration chosen, sequence-design is a process of finding a sequence that gives maximum depth 
to the funnel. It is worthwhile investigating whether a mutually consistent treatment can be given to 
these two processes - namely, (i) identification of a target configuration which is not very sensitive to 
changes in the contact-energy scheme and (ii) designing a sequence that ensures maximum energy- 
gap between the target configuration and the others. 

In this paper, we describe an exact enumeration-based search algorithm in the combined sequence- 
structure space for identifying the target (native) SAW configuration as well as its designing se- 
quence on a diamond as well as on a square lattice. The motivation is to evolve an algorithmic 
criterion for identifying the best energy scheme from among the set of schemes being considered. 
We show that the target configuration (a 25-step compact SAW on a diamond lattice llBIl and a 
24-step compact SAW on a square lattice) is not sensitive to the different energy schemes consid- 
ered, whereas the designing sequences are. The native configuration on a diamond lattice is seen 
to consist of two 'sheets' (Fig[T]) indicating the appearance of secondary structural elements even 
for such short chains. From the Density of States and the numerical values obtained for the zeros 
of the partition function, we calculate the temperature-dependent specific heat for all the designing 
sequences and the energy schemes. However, ground-state degeneracy (i.e., the number of target 
configurations) increases for shorter walks. We have also that the target configurations of shorter 
(23-, 21- and 17-step) walks on a diamond lattice show sheet-hke secondary structures. 



Methods 

A. Self Avoiding Walks (SAW) on a lattice Hi. 

An ordered set of points on a given lattice of coordination number z, called a 'chain' or a 'walk', 
can be specified either by listing the coordinates of all the points or by specifying the starting point 
and the sequence of direction-vectors between consecutive points of the walk. It is a SAW if the 
points in the list are all distinct, which is to say that the walk does not visit a site already visited 
once. 

Different lattices have different sets of direction-vectors. For example, a cubic lattice has the set 
of six direction-vectors {(±1, 0, 0), (0, ±1, 0), (0, 0, ±1)} whereas a diamond lattice has two sets 
of direction-vectors, namely ±{(1, 1, 1), (—1, —1, 1), (1, —1, —1), (—1, 1, —1)} , which have to be 
used alternately in going from site to site. 

SAWs are known to have certain universal properties, independent of the lattice in which they are 
embedded. For example, the total number of iV-step SAWs on any given lattice has an asymptotic 
scaUng form Ist], 

Zn = awiV" (1) 
where 7 is a universal exponent whose value depends only the dimensionality of the lattice, /i is 



the effective coordination number 111411 and is a lattice-dependent proportionality constant.The 
average 'size' of an A^-step SAW is specified either in terms of the mean squared Radius of Gyration, 
{Rg{N)), or equivalently in terms of the mean squared end-to-end distance, {Rl{N)), which are 
again of the form, 

{Rg{N)) =bNN^'' (2) 
{RliN)) = cnN^" (3) 



where 1/ is a universal exponent 111411 and, and cn are lattice dependent proportionality constants. 
In the limit of large N, the ratio, b^/cN, is known to converge to the values, 0.14 and 0.16, in two 
and three dimensions respectively. 
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An A^-step SAW configuration may have sites that are nearest neighbors on the lattice but far away 
along the chain; such 'non-bonded' pairs of nearest neighbors are referred to as contacts. If Z^r ^ 
denotes the number of A^-step SAWs having m contacts, we can also write 

AI 

Zn — ^ Zjv,™ (4) 

where AI is the maximum number of contacts an A^-step walk can have. Labeling the sites of a 
configuration in a sequential order leads to a concise list of contacts - namely, {0, fc) | j,k = 
1, 2, 3, A^; J < fc} - called its co«facf-ma/5. 



B. HP-model 0]. 



In this model, conformations of proteins are assumed to be represented by lattice SAWs whose 
sites are occupied either by hydrophobic (H) or by polar (P) residues. So, an HP-protein is a point, 
{Ca/'j (?}, in the combined sequence-structure space where Cj^ is the backbone SAW configuration 
to which is assigned the binary sequence a = {aj \ aj = H or P; j = 1,2,3,..., A^}. We treat an 
amino acid residue (H or P) in the chain as a structureless bead and refer to it simply as a 'monomer' 
of type H or P. 

Each monomer is connected (or bonded) to two other monomers in the chain and hence can make 
a maximum of {z — 2) contacts; however, at either end of the chain, it can have a maximum of 
(z — 1) contacts. Interactions between monomers contribute to the energy, E{Cn,<j), of an HP- 
protein. Bond-energy contribution to E{Cn, f?) can be ignored as it is a constant independent of 
the chain-configurations; so, it is only the interactions between non-bonded nearest neighbors in the 
chain that contribute to E{Cn, f?)- Different configurations can have different sets of contacts and 
hence can have different values of E{Cn, <?)■ In terms of a basic unit of energy, cq, we may write 

E{CN,a)^ CjkUia,,ak) (5) 

(i,fc>j+i) 

where Cjk is the contact-map given by 

„ f 1 for non-bonded nearest neighbors i and j , 

^J*^ \ othewise 

and U{(jj,(Tk) denotes the energy associated with the contact (j, k) whose value depends on the 
type of monomers (H or P) in contact - namely, U {aj,ak) = {uhh, UHP,upp}. In other words, 
U{aj,ak) depends not only on the configuration through its contact-map but also on the HP- 
sequence considered. 



Based on a statistical analysis 11511 of the various interactions between the amino acid residues 



of real proteins, the values of U{<Tj, ak) have been shown lll6ll to satisfy the inequalities, uhh < 
uhp < upp and {uhh + upp) < 2uhp- Together, these inequalities imply the dominant role 
of hydrophobic interaction in getting the H-type of monomers together and also a general tendency 
of the two types of monomers to segregate. The number of two-lettered HP sequences leading to a 
unique ground state, for a given scheme of energy values, is taken to define the 'designability' of a 
model. 

The simplest energy scheme is the one in which the only interaction considered is that between 
the H-type monomers, i.e., {uhh, uhp, upp) = (—1, 0, 0) IitIi . It is equivalent to saying that the 
minimum energy configurations, referred to as the ground states, are those with maximum number 
of contacts. This model brings out the gross qualitative features that are comparable to those of 
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real proteins; however, it leads to highly degenerate ground states with small number of 'designing' 
sequences (i.e., sequences with unique ground state). 

A modified version of this simple scheme is called the shifted HP-model in which the it's are 
assigned the values, {uhh, uhp, upp) = (—2, 1,1) llSll . Yet another widely studied scheme is 
defined by the set of values, {uhh ,uhp, upp) = (—2.3, —1, 0) ifTisll . It is not necessary that the 
ground states in these schemes should correspond to maximally compact {i.e., maximum number of 
contacts) configurations. However, they lead to less degenerate and better designable structures. As 
pointed out by Pande et al. ifioll . similar interactions lead to similar degrees of designability while 
different interactions lead to different patterns. 



C. Target structure. 

By definition, a target structure is a SAW configuration, usually a maximally compact one, that 
is assigned minimum energy by a set of HP-sequences within a contact energy scheme, say U{= 
{uhh, Uhp, Upp}). The value of the minimum energy realized depends on the configuration, and 
is given by 

N\ 

Eo{CN,nH;U)=mm{E{CN,aj)\j = 1,2,3,.. .,M}; M = — (7) 

nH\[N ~ nny. 

where tih is the number of H-type monomers that are required to be part of the HP-chain. Let 
b{CN, TiH', U) denote the the number of HP-sequences assigning minimum energy, to the configura- 
tion Cn- Then, it is more designable than another configuration, say C^, having the same minimum 
energy if b{CN,nH', U) > b{C'j^,nH; U). Given a set of configurations, {Cn}, those having the 
energy Eg = min{i?o} are the ground state configurations each of which in turn is characterized by 
a certain value of h. Configurations of interest are those, {C^}, that are characterized by the values 
{Eg, bg) where bg = max{&}. 

It is not necessary that the set of designing sequences be the same for all the configurations, C^. 
Also, the sets of designing sequences as well as the set of designing configurations, {C^}, may be 
different for different contact energy schemes. In order to confirm this explicitly, we have considered 
the following contact energy schemes: 

Ui = (-1,0,0) (8) 
U2 = (-7,-3,0) (9) 
t/3 - (-7,-3,-1) (10) 
Ua = (-7,-3,1) (11) 



Ui is the standard energy scheme mentioned earlier 111711 ; [/2 is just three times the scheme 
(—2.3, —1, 0) 11161 , integerized; includes an attractive PP-contact energy, and is not related to 
U2 by any numerical operation such as integer addition or multiplication. All these three schemes 
satisfy the inequalities ||l6l] mentioned in the last section. The scheme has a repulsive PP-contact 
energy, and violates one of the inequalities namely, [uhh + upp) < 2uhp', this is done just to see 
what difference it makes vis-a-vis the other schemes. 



D. Stability. 



Levinthal's paradox lfl9ll refers to the observation that a denatured protein folds spontaneously to 
its native state on very short time scales even though a sequential search for its native state would 
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take forever to complete. It was then thought that folding proceeds along specific pathways through 
intermediate structures. But, Anfinsen's yj experiment demonstrated that folding is independent of 
pathways, suggesting thereby a funnel-like free energy landscape in which the native state lies 
at the bottom. While the depth of the funnel corresponds to the energy-stabilization of the native 
state with respect to a denatured state, the width of the funnel corresponds to the conformational 
entropy of the chain. This picture provides a conceptual framework for a fast folding kinetics. 
The existence of a large energy gap between the native state and any other conformational state is, 
therefore, considered as a prerequisite for a model protein . 

As described in the previous section, we may have a set of ground state configurations, {C^}, 
all having the minimum energy Eg and designable to the same degree; that is to say, the number of 
HP-sequences associated with each configuration is the same, namely bg. The sequences themselves 
may be different for different configurations. 



1. Non-degenerate ground state: 



If there is only one ground state configuration, C^, and the associated set of HP-sequences (all 
giving the same energy, by definition, to C^), it is necessary to check whether it is at the bottom of 
a funnel in the energy landscape. 

Since Cf^ is the unique ground state, all the other configurations are excited states whose energies 
will depend upon the HP-sequence fitted to them. For every sequence G S^, an average excited 
state with energy {Eg{sg)) and the corresponding standard deviation Ag{sg) can be obtained. The 
gap, {{Eg{sg)) — Eg) is a measure of the funnel-depth and Ag{sg) is a measure of the funnel- width. 



The ratio of these quantities 112 III 



fn,sg) ^ ((^.^(^g)) ^^> . f^c%) ^ max{/(C^,.,)} (12) 



can be taken to be a single parameter characterizing the funnel. Sequences having the maximum 
gap-parameter value, /, may be deemed the best designing sequences. 



2. Degenerate ground state: 



If there are more than one ground state, then a priori all are equally probable and so the funnel is 
likely to have so many wells at the bottom (a flat bottom, in the coarse-grained sense). With each of 
the configuration, C^, we have a set of sequences S^{Cf^) all giving it the minimum energy Eg. If 
there is no sequence common to all the sequence-sets (i.e., flcs^ S^{Cf^) = (f>), then we may identify 
those configurations for which the gap-parameter, /(= max{/(C^)}), is the maximum. On the 
other hand, if 0^9^ S^{Cf^) ^ 0, then this procedure cannot reduce the degeneracy of the ground 
state; all the degenerate configurations have the same set of common designing sequences. 

It is possible that different contact energy schemes give rise to different sets of ground states. 
Since the estimate of a contact energy is statistical in nature llisll . it is not unreasonable to speculate 
that there may be ground state configurations that are insensitive to changes in the energy scheme. If 
such is the case found by repeating the above procedure for different energy schemes (Eq.dSlfTTt. for 
example), then degeneracy is reduced even further In other words, we try to identify those configura- 
tions for which the gap-parameter, /(= max{/(C^)}), has a set of maximum values corresponding 
to the set of contact-energy schemes. 

It is also possible that there are configurations designed by sequences that are mutually in- 
verse (mirror-image) to each other; for example, HP/i{H2P)2HPH^P2HP-iH2 is inverse to 
H2P7iHP2H'iPH{PH2)2PiH . In that case, we may check whether the configurations can be 
mapped into one another by a simple permutation and relabeling of the direction-vectors; if so, 
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they are equivalent configurations. If this procedure leads to just one configuration, then we have a 
unique ground state along with its set of designing sequences; otherwise, this procedure reduces the 
ground-state degeneracy to a maximum extent. 

3. Specific Heat: 

Given a sequence and an energy scheme, we use the enumeration algorithm to have an exact count 
of the number of configurations having an energy, E; The integerized forms of the energy schemes 
chosen, (Eq.dSlfTTTi). ensure that E is an integer between zero and the minimum value —Eg. The 
total number of an N-step walks, also called the microcanonical Partition Function, is then given by 

Zn{U) = Y. Zn.e{U) (13) 

£=0 



which is the same as Eq.(|4|i for U ~ Ui l24tl . The numbers, {Z ^ ,e{U)} , are also known as the 
Density of States (DoS). At finite temperatures, T, the DoS will be weighted by the Boltzmann 
factor, e^l^\ since E is an integer, the weighted sum, Zn{U,w), called the canonical Partition 
function is actually a polynomial of degree Eg in the variable w{= e^/^). We can write 

ZN{U,w)^Y[iw^WE) (14) 

E=0 

where wes are the (complex) roots of the polynomial equation, Zn{U, w) — 0. These roots will 
occur in the form of complex conjugate pairs because Zn{U, w) is a real number The specific heat 
at any given temperature is then given by 

C'(U,w) = w[\n[w)\ — Iw — 1 (15) 

In terms of the real, u, and the imaginary, v, parts of w, we have 

„ \(2W — — U^(W — U,-)^l 1 /rr 

ciu,w)^wMw)rJ2- — u ^2!. 212 ; ^ = «/ (16) 



Fast folding sequences have been observed II25I1 to give rise to two different peaks in the specific 
heat of which the one at lower temperature is expected to be sharp and prominent. 

It must be mentioned here that the specific heat, C{U, T), can also be computed directly from the 
DoS using the expression, 

^ TrE=oZ{N,E)^E/T 

and it will be the same as C(C/, w). However, partition function zeros provide additional corroborat- 
ing evidence for considering a peak in the specific heat as indicative of a phase transition, especially 
when we study small systems. 



Results 



1. Algorithm: 
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Using an exact enumeration module for generating SAWs on any given lattice, we first dump all 
N-step SAWs having maximum number of contacts. For a given energy scheme, U, and a binary 
sequence having nn ones, it is the set of contacts that determines the energy of a walk; so, we 
identify and collect all the walks with distinct contact-maps. 

In order to identify the target configurations {i.e., walks, to which maximum number of HP- 
sequences assign minimum energy), we generate binary sequences having a specified number of 
ones ('H's) and find the energy they assign to the maximally compact walks, taken one by one, 
by using their contact-maps. This way, we have a count of binary sequences that assign minimum 
energy attainable for each walk and a list of those walks having global minimum energy; from this 
list, we choose walks with which maximum number of (designing) sequences are associated. These 
are the (designable) target configurations for a given energy scheme, say {C^v.nff (g)}!7- 

By repeating this procedure for all the specified energy schemes, we have sets of designable 
configurations and their associated sets of designing sequences. If there are configurations common 
to all these sets, say {CN.nnig)} = ^u{CN,nH{9)}u, we consider only them for further analysis; 
otherwise, we have to consider all these sets one by one. As described in the last section, we use the 
gap-parameter,/, for reducing the degeneracy of the ground state as far as possible: 

From the set of sequences designing a given configuration, say Cj^''„ (g) £ {Cjv.njf (5)}, we 
identify the ones with maximum values of f{U) for the energy-schemes; it is possible that the same 
sequence has maximum value of for more than one energy-scheme. We do it for every configuration 
in the set; the set of sequences and their corresponding / values may be different for different 
configurations. For example, if there are K configurations under consideration (j — 1, 2, • • • K) 
and four energy-schemes, then for each j, we have four sequences (same or different) and their 
corresponding / values, say f'^^'> = {fi\f2\f3\fi^}- We choose those configurations and 
their associated sequences to which correspond the maximum of the set {/'^'}. If there are more 
than one such configuration, then we have to decide whether they are equivalent with respect to a 
permutation and relabeling of the direction-vectors. This way, the ground state degeneracy can be 
reduced as much as possible. 



2. Diamond Lattice walks: 



In the case of diamond lattice, we have two sets of direction-vectors, namely {«/(.?) | j = 
1,2,3,4} = {(1,1,1), (-1,-1,1), (1,-1,-1), (-1,1,-1)} and {«„(j)} = {-«/(i)} which 
have to be used alternately in going from site to site — i.e., if the A'*'' step is taken in the direc- 
tion vi{j), then the next step will have to be taken in one of the three directions vjj{l), I 7^ j after 
ensuring that this step avoids previously visited sites. A recursive procedure for adding a step in a 
systematic manner can enumerate all SAWs up to a specified maximum length. 

Also, any given configuration can be coded in terms of the direction-labels such as, for example, 
1234231242313123243212132 which codes a 25-step walk. Using such a code and assuming that 
the walk starts from position ro(from the origin, for example), the positions of the subsequent sites 
in the walk can be generated by knowing the direction-type (I or II) of the first step. We have fixed 
the first two steps in the directions u/(l) and u//(2) respectively. 

Using such exact enumeration procedure, we have generated 25-step SAWs on a diamond lattice 
II22II and found that there are 270 maximally compact walks with 12 contacts each, of which 135 
walks have distinct contact-maps. In TABLE IH we have presented the number of walks 1231] and 
their average end-to-end distance as a function of the number of contacts. 

By fitting all possible 26-bit binary sequences (13 ones and 1 3 zeros) to each of the 135 maximally 
compact and distinct configurations, we shortlisted designable configurations for each of the energy 
schemes given in (Eq.dSlfTTTi) and presented the data in TABLE For a given energy scheme, the 
number of sequences associated with each of the 40 ground state configurations is the same, and 
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m 


^25 (m) 


(r25(m)> 





18409190935 


15.3801 


1 


16564041012 


13.4874 


2 


9902194362 


12.0026 


3 


4583935502 


10.5940 


4 


1791649410 


9.3488 


5 


606531194 


8.2332 


6 


180088606 


7.2606 


7 


44257966 


6.4466 


8 


8323860 


5.7982 


9 


1068592 


5.2175 


10 


67036 


4.9458 


11 


872 


5.5048 


12 


270 


3.7446 



TABLE I: ^25(01): Number (xl2) of 25-step walks making m contacts; {r25(m)): average end-to-end dis- 
tance. 



u 




E, 


So 




40 


-8 


630 


U2 


40 


-59 


420 


U3 


40 


-62 


420 




40 


-56 


31824 



TABLE II: U : Energy scheme (Eq. (l8ll ft : Z!^^ [U): Number of 25-step ground state configurations correspond- 
ing lo U; Eg] Ground state energy; Sg: Number of designing sequences. 



u 


Sequence 


f 


Ui 


HzP.2HP2H2P2HP2H2PiH2{PH)2 


4.2126 


U2 


HiPH2P3H2PHP2{HP)2PHP{PH)2 


3.9946 


Ui 


{HP)2H2PiH2PiH-iP2H2{PH)2 


4.1622 


Ui 


H3PH2P3H2PHP2{HP)2PHP{PH)2 


3.7132 



TABLE III; f; maximum gap-parameter. Subscript of H or P denotes the number of times the corresponding 
symbol is repeated. 

is the maximum. The sets of ground states, though 40 in number for all the schemes, are different 
except for one configuration, 1231212313124231312323123. This is the configuration which is 
also the ground state for all the four energy schemes and hence is the most designable ground state 
configuration. As shown in Figd] this configuration is clearly seen to consist of two sheet-like 
secondary structures. 

Among the sets of sequences designing this configuration, for the energy schemes considered, 
We found 420 sequences common to these sets. In other words, these 420 sequences design the 
configuration irrespective of the energy scheme considered and hence may be referred to as the best 
designing sequences. We used these sequences and the original set of 135 distinct configurations so 
as to address the thermodynamic stability of the unique ground state configuration obtained (FigdJ. 

As described above, for every energy scheme (C/i to Ua), we find the gap-parameter associated 
with each sequence by identifying which of the 135 configurations correspond to 'excited' states. 
From the list of the gap-parameters associated with the corresponding sequences, we can choose the 
sequence with the maximum gap-parameter. We present the best sequences obtained in TABLE Hill 
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FIG. 1: Two different views of the most designable configuration, 1231212313124231312323123, that is com- 
mon to all the four energy schemes considered (Eg. (181 111 ) . As is evident from view (b), the structure consists 
of two 'sheets'. 



We see that the sequences corresponding to U2 and U4 are the same and so, we have three distinct 
sequences. This result suggests that the basic energy scheme also plays a role in ensuring a large gap 
between the ground state and the excited state configurations. To further confirm that the schemes 
U2 and U4 lead to the same designing sequence, we have rerun the exact enumeration algorithm for 
computing the Density of States (DoS) for these sequences and the energy schemes. 

We have computed specific heat as a function of temperature for the three distinct sequences, say 
Si , S2 and 5*3, and the energy schemes Ui — U4; the results are shown in Fig|2] The insets show the 
complex zeros, we, of the partition function defined in Eq.(fT4li. 

As can be seen from the figure, all the three sequences show a sharp peak at T ~ 0.4 for the 
energy scheme Ui. They also show more or less the same behavior for the energy schemes U2 
and U4 implying thereby that these two schemes, except for a scale-shift, are probably the same. 
Sequence Si has a hint of a shoulder at T ^ 1.0 (Fig|2a)). On the other hand, sequence 6*2 shows a 
prominent peak at T ^ 1.5 for the energy scheme U3 (Fig|2jb)) whereas only a shoulder at the same 
value of T for U2 and U4. Sequence S3 shows a shoulder at T 1.4 only for the scheme U3. 

The broad peak seen at T ~ 2.3 is indicative of a coil-to-globule transition whereas the ones 
observed at T ~ 1.0 — 1.5 are indicative of a globule-to-native folding transition. According to 



Socci and Onuchic II25I1 . fast folding sequences give rise to two different peaks in the specific heat 
of which the one at lower temperature is sharp and prominent. This observation implies that the 
sequence S2 with the energy scheme U3 is a fast folder as compared to the others. 

It must be mentioned here that configurations with less number of contacts for these sequences and 
the energy schemes may have lower energy than the ground state values reported in TABLE II. We 
summarize the results in TABLE HV] From these tables, it is clear that the designability {i.e., number 
of designing sequences) of these lowest energy, less compact (only 10 contacts) configurations is 
much less than that of the most compact (12 contacts) configurations. 



s 




E,{U2) 


EAU3) 




CC 


51 [2] 


-9(4) 


-63(4) 


-64(4) 


-62(66) 


1231212313124231432421321 


5'2[13],S3[25] 


-9(33) 


-66(9) 


-66(9) 


-66(9) 


1231 24234342 1 23 1 2 1 243242 1 



TABLE IV: Diamond lattice walks. Number of contacts = 10. Numbers within square brackets denote the 
number of distinct configurations; those within curved brackets denote the number of designing sequences. 
Eg{U): minimum energy for the scheme U ; CC: common configuration 



A question arises whether shorter chains also have a unique target configuration such as the one 
shown in FiglTJ so, we have repeated these calculations for walks consisting of 23, 21 and 17 steps 
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N 


M(N) 




25 


12 


1 


23 


10 


2 


21 


9 


15 


17 


6 


68 



TABLE V: Zf^: the number of maximally compact N-step minimum-energy configurations on a diamond lat- 
tice, which are common to all he energy-schemes, Ui — U4. M(N) is the maximum number of contacts for 
walks of length N. 

and have presented the data in TABLE [V] We note that the number of minimum-energy configura- 
tions common to all the energy-schemes Ui — U4 (i.e., the ground state degeneracy) increases for 
shorter waUcs. 

Ti-step walk: There are two 23-step target configurations (TABLE IVH for which we carried 
out the sequence-design process and the results are presented in TABLES IVIII - IVIIII There are 
two distinct sequences designing the configurations, and Cfg^'' - one for the energy-schemes 
Ui & C/3 and the other for the energy-schemes U2 & Uj^. It is interesting to note that the sequences 
that design Cfg^' are the inverse of those designing Cfg^'; for example, the sequence Si designing 
Cff forthe energy-schemes [/i & [/a is the inverse {i.e., mirror-image) of the sequence designing 
Cfg^'' for the same energy-schemes. Even their respective gap-parameters have roughly the same 
value. Moreover, the configurations Cfg^' and Cjs^'' shown in Fig[3]can be mapped into each other 
as follows. 

The reversed configuration code of Cfg^' is 23123231214213132121321; a clockwise permuta- 
tion of the direction-vector labels (123) — )• (231) and relabeling them as (123) leads to the code 
12312123134132321313213 which is exactly that of Cjs^'. Hence, we have effectively a unique 
target structure. 

2\-step walk: There are fifteen 21 -step configurations common to all the energy-schemes for 
which we carried out the sequence-design process and the results are presented in TABLE |IX] 
We note that the same squence, H2PH^PH2PaH2P3{HP)2, designs the five configurations C|p\ 
^21^^ ^21^' ^21^ ^21^ for the energy-schemes U2 — U4; the gap-parameter / is also the same. 
The sequence designing Cff ^ for the energy-scheme Ui has the highest value of / (= 6.4846). 
Among the fifteen structures, Cj^^"* has the highest set of / values and hence may be chosen as the 
target structure (FigH). 

ll-step walk: There are sixty eight 17-step configurations common to all the energy-schemes for 
which we carried out the above-described sequence-design process. We found two configurations 
being designed by the same sequence {PHP^H^{PH)2{HP)^, for all the energy-schemes) with 
maximum values for the gap parameter - data presented in TABLE [X] Further, we note that these 
configurations are mappable into each other by the the relabeling of the direction vectors 3 < — > 4 
and hence may be considered equivalent (Fig|5]l. 

To sum up, in the case when we have ground-state degeneracy, we can choose a target structure 
by adopting either of the two methods or both - namely, (i) checking whether there are configura- 
tions whose designing sequences are inverse of each other and whether these configurations can be 
mapped into each other by a permutation and relabeling of the direction- vectors, and (ii) choose the 
configurations that have highest set of / values. This way, the ground-state degeneracy can be mini- 
mized implying thereby an efficient approach towards realizing a funnel-shaped energy-landscape. 

3. Square Lattice walks: 

Just to check whether the algorithm described above for diamond lattice walks also works for a 
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zuu) 
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Ui 


34 


-8 


10 


U2 


34 


-56 


10 


Us 


34 


-58 


10 




2 


-54 


1001 



TABLE VI: Data for 23-step walks on a diamond lattice; definitions are the same as the ones in TABLE HH 
The two configurations common to all the energy-schemes are, (i)C^^^': 12312123134132321313213 and (ii) 
C?P^: 12312123131241213232132. 



u 




Sequence 


f 


Ui 


Si : 


HPi{H2P)2HPHsP2HPsH2 


4.1783 


U2 


& 


: HPHP2H2P2HPH{P2H2)3 


3.5876 


Us 


Si : 


HPi{H2P)2HPH2.P2HPiH2 


3.8694 


Ui 


S2 


: HPHP2H2P2HPH{P2H2)3 


3.2946 



TABLE VII: Cl^^: 12312123134132321313213; definitions are the same as the ones in TABLE|lIIl 



U Sequence f 

f/i Si : H2PiHP2H'iPH{PH2)2PiH 4.1776 

U2 S2 ■■ {H2P2)3HPHP2H2P2HPH 3.5860 

Ui Si : H2P3HP2H3PH{PH2)2PiH 3.8759 

Ui S2 ■■ {H2P2hHPHP2H2P2HPH 3.2959 

TABLE VIII: C£f ^ : 12312123131241213232132; definitions are the same as the ones in TABLE|lIll 



square lattice (same coordination number), we have considered 24-step walks on a square lattice. 

There are a total of 1081 24-step maximally compact (16 contacts) SAWs on a square lattice with 
distinct contact-maps. By fitting all possible 25-bit binary sequences (12 ones and 13 zeros) to each 
of these configurations, we found one configuration, 112234323341434412141123 I26I1 . that is the 
minimum energy configuration for all the four energy schemes Ui — U/^. 

For all these energy schemes, we found just 4 sequences common to the sets of sequences de- 
signing this configuration. We estimated the gap-parameters for these sequences by using the whole 
set of maximally compact walks, and found that just one of these four sequences has the maximum 
gap-parameter, /, for Ui — U4 (TABLE lXib . The corresponding ground state configuration is shown 
in Figlgta). 

Using the DoS, computed for this designed sequence and energy schemes, we have calculated the 
specific heat and presented the data in Fig|6jb). The prominent first peak seen for U3 implies that 
this designed sequence is a fast folder only for this energy scheme. 

There are six distinct less compact configurations with 14 contacts having less energy (= —78) for 
this designed sequence than the ground state energy (= —76) of the maximally compact walks, but 
only for the energy scheme U4. While the number of sequences designing the maximally compact 
ground state configuration (Fig|6la)) is 220, it is much less (= 32) for these less compact walks. In 
other words, the latter are not only less designable for U4 but are actually excited states for the other 
energy schemes. 

We also found that the designability also depends on the number of ones (H type) in the sequence. 
For example, ground state configurations for a designed sequence consisting of 13 ones are shown in 
Fig|7l;a,b). Both are designed by the same sequence HP2HPH2P2HPH2P2HPH2P2HPH2. The 
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ground state displayed in Fig|7Ja) is only for the scheme U4 with a designability of 195 sequences, 
whereas the one displayed in Fig|7jb) is common to all the energy schemes. But the designability of 
the latter is just one sequence for the schemes Ui — U3. We also see that these configurations don't 
have a clear hydrophobic interior as compared to the one, Fig|6ja), designed with 12 H-types. 

Conclusions 

We have presented an exact enumeration study of how the sequence-design of maximally com- 
pact lattice SAWs is influenced by the various contact-energy schemes chosen; the motivation is 
to evolve an algorithmic criterion for identifying the best scheme from among the set of schemes 
being considered. Moreover, since the various contact-energy schemes considered in the literature 
are based on statistical estimates of the average effective inter-residue interactions of real proteins, 
it is of interest to know if the target configuration and its designed sequence is stable with respect to 
changes in the contact-energy scheme. In fact, the many-to-one mapping of amino acid sequences 
onto the set of all native structures implies that a given native structure should not be very sensitive 
to the set of contact-energy schemes underlying the set of sequences designing it. 

We have considered four different schemes, Ui — U4 (Eq.dSlfTTTl. of which the last one has been 
deliberately chosen not to satisfy the inequality {uhh + upp) < 2uhp that the first three satisfy. 
The enumeration data for maximally compact 25-step walks on a diamond lattice, as well as for 
24-step walks on a square lattice, suggest the following. 

(a) The target structure {i.e., the most designable, maximally compact configuration) is the same 
for all the energy schemes considered: (i) On a diamond lattice, such a configuration consists of two 
sheet-like secondary structures as is evident from FigH] (ii) On a square lattice also, it is independent 
of the energy scheme and is shown in Fig|6ja), which is better designable with 12 rather than with 
13 H-type monomers Fig|7Ja). /.e., designability of the target structure, rather than the bare structure 
itself, could very well depend on the proportion of H-type monomers in the sequence. 

(b) The designing sequence may depend on the contact-energy scheme: (i) On a diamond lattice, 
we find that the designed sequences for U2 and U4 are the same, implying thereby that these two 
energy schemes are the same from sequence-design point of view, whereas those for Ui and U3 are 
different (TABLE III); (ii) On a square lattice, however, the designed sequence (12 H-type) is the 
same for all the energy schemes (TABLE V). 

(c) We find that the maximally compact ground state configuration is better designable than a few 
less compact configurations that are assigned lower energy by the designed sequence. 

(d) On the basis of an exact count of the Density of States, we have calculated the temperature- 
dependent specific heat for these designed sequences and the four energy schemes, (i) In the case 
of diamond lattice walks, sequence S2 with scheme U3 seems to be the best fast-folder because it 
has a prominent low-temperature peak in the specific heat that corresponds to the compact-to-native 
transition; (ii) In the case of square lattice walks also, the designed sequence (common to all energy 
schemes) with U3 seems to be the best folder for the same reason. These specific heat data imply 
that the energy scheme U3 is preferable to the other three. 

(e) We have confirmed that this algorithm works for shorter walks as well by presenting data for 
23-, 21- and 17-step walks on a diamond lattice. We have been able to identify a unique target config- 
uration by taking care of the ground-state degeneracy as described in the last section. Interestingly, 
all these shorter target configurations also show sheet-like secondary structures. 

To sum up, we have presented an algorithm for identifying the target configuration that is the 
minimum energy configuration for a given set of contact-energy schemes; the process is implicitly 
sequence-dependent. The best designing HP-sequence for this configuration is then the one for 
which the gap-parameter has the maximum value. That an implicit sequence-dependence of a folding 
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process ensures a realization of a funnel-like energy landscape has been demonstrated recently by 



Chickenji et al 112711 . In particular, they have shown that short fragments excised from proteins 
suggest sequence-dependent local restrictions in the conformational space that are also relatively 
free from the inaccuracy of the local interactions. While their method assembles appropriate short 
fragments into a target configuration, our method is a sequence-dependent search procedure for a 
target configuration that is not strongly dependent on the contact-energy scheme. There is a need for 
a more detailed study of our method vis-a-vis others. 

SLN thanks Artur Baumgaertner for helpful discussions and a critical reading of the manuscript. 
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FIG. 2: Specific Heat as a function of temperature. Insets show the zeros of the partition function Eq.([T4) for 
the energy schemes U2, U-j, and Ui. The sharp peak at T ~ 0.4 (black) corresponds to Ui. (a) Sequence Si: A 
Shoulder is seen for energy schemes U2 (red) and U4, (green) at T ~ 1.0. (b) Sequence ^2: A shoulder is seen 
for energy schemes U2 (red) and U4, (green) at T '-^ 1.5; a prominent peak at about the same value of T is seen 
for energy scheme f/3 (blue), (c) Sequence S^: A shoulder is seen at T ~ 1.4 only for the scheme f/3 (blue). 
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FIG. 3: 23-step configurations Cls^' and Cfg'^^ that are common to all the four energy schemes considered 
(Eq. dSIl it ) . As is evident from these figures, they consist consists of two 'sheets' . 




FIG. 4: 21 -step configuration CjP' that is conmion to all the four energy schemes considered (Eg. JSIl ll ) and 
has the highest set of / values. Two 'sheet' -like structures can be seen. 




FIG. 5: 17-step configuration Cfj that is common to all the four energy schemes considered (Eg. JSIl ll ) and 
has the highest set of / values. Again, 'sheet'-like structures can be seen. 
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C^{^': 123413432121341312313 
U Sequence f 

C/i H2PAHP2H2P2{HP)2HiP 5.1075 
U2 {PH)2HPHP2H2PaHP2Ha 4.3989 
U'i {PH)2HPHP2H2PiHP2Hi 4.7583 
U4. PH2P2HzPH2PiHP2H'iP 3.8996 

C^f^ 123423132131431212343 
Ui H2PiH3P{PH)2P2H3PHP 4.7633 
f/2 H2PiHiP2HPH2PH2P2HP 4.2390 
(73 H2P4H3P2HPH2PH2P2HP 4.9193 
HiP2HP3H{HP)2H2P2HP2 4.7140 
C^f^ 123124121343212412312 
;7i {PH)2PiH2PHP2HzP2H-iP 6.1679 
f/2 H2PHzPH2PiH2Pi{HP)2 6.0343 
[/3 H2PHzPH2PiH2Pi{HP)2 6.5275 
f/4 H2PHzPH2PiH2Pi{HP)2 5.2154 

Cff : 123124121343212432132 
[/i H3P2H3P2H2P2HP3H2P2 6.3509 
[/a H2PHzPH2PiH2P3{HP)2 6.0343 
[/3 H2PH3PH2PAH2Pa{HP)2 6.5275 
f/4 H2PH3PH2PiH2Pi{HP)2 5.2154 

C^f^ 123124131242313413213 
T^i P2H2P3HP2H2iP2Ha)2 6.4846 
f72 H2PH3PH2PaH2P3{HP)2 6.0343 
[/a H2PH3PH2PiH2P3{HP)2 6.5275 
f/4 H2PH3PH2PiH2P3{HP)2 5.2154 

C^f^: 123124131242313423123 
f/i {PH)2P3H2PHP2H3P2H3P 6.1679 
f/2 H2PH3PH2PiH2P3{HP)2 6.0343 
f/3 H2PH3PH2PiH2P3{HP)2 6.5275 
f/4 H2PH3PH2PiH2P3{HP)2 5.2154 
C^f^ 123243132131231412323 

C/i P2HP3H2PH3P2H2P2H3P 4.2178 
f/2 P2HP3H2PH3P2H2P2H3P 4.1126 
f/3 P2HP3H2PH3P2H2P2H3P 4.3714 
f/4 HiP2HP3H2{PH)2{HP2)2 4.7140 
Cff ; 123213232141321312313 
f/i {HP)2PHPH2P2{HP)2H3P2HP 5.6569 
C/2 {HP)2PHPH2P2(HP)2H3P2HP 4.7329 
C/3 {HP)2PHPH2P2{HP)2H3P2HP 5.1549 
f/4 HPH3P2H2P3{PH)3P2H2 3.8735 



C^f^ 124123121434212342142 
f/ Sequence f 

f/l H3P2H{H2P2)2HP3H2P2 6.3509 

f/2 H2PH3PH2PiH2P3{HP)2 6.0343 
f/3 H2PH3PH2PaH2P3{HP)2 6.5275 
f/4 H2PH3PH2PiH2P3{HP)2 5.2154 
C^f ^ 121234314321241312313 
f/i HPH3P3HP{PH)3P2HPH2 5.2511 

f/2 P2H2P{PH)2P2H2P2HP2H4. 4.0692 

f/3 HPH3P3HP{PH)3P2HPH2 4.9029 
f/4 PH3P2HP{PH)2P3H2PH3P 3.7262 

C^{"^ 121343413431231413414 
f/i H2P2H2P4HP2HPH3{HP)2 4.4721 
t/2 H2PH3P4HP2HPH3P2HP 3.9009 
t/3 H2PH3P4HP2HPH3P2HP 4.2762 
f/4 P2HP2H{HP)2H2P3HP2H4 4.0415 
Cff^^ 121321234231321313242 
f/i H2Pi{HP)2H2PH{HP)2{PH)2 3.6181 
t/2 PHP2H5P2HP2{P2H2)2 3.6856 
t/3 H{HP)2H2P3H2PHP2HP2H2P 4.6512 
t/4 HiP2HP3H2{PH)2{HP2)2 4.7140 

C^f^^ 121321231242313123132 
f/i PHP2H3{PH)2P2H{HP)2{PH)2 4.7329 
f/2 PHP2HT,PHP2H2Pi{PH)2 4.0565 
f/3 PHP2H5PHP2H2P3{PH)2 4.4721 
f/4 PHP2H3PH3P2H2P3{PH)2 3.5820 

C^f''^ 121321241232314124132 
t/i PH3P2HP4HSP3H2P 4.5210 

f/2 PH3{P2H)2PH2PH2P3H2P 4.2083 
f/3 PH3P2HP4H5P3H2P 4.4488 

f/4 H2PH{P2H)2PH2PHP2HPH2P 4.3894 

c^J^'^^ numA'ii'iuiuui'ii 

Ui H4P2HP{P2H)2{H2P2)2 4.5518 
f/2 PH3P2HPH3P2{P2H)2HPH 4.4783 

f/3 H4P2HP{P2H)2{H2P2)2 5.2152 

f/4 H3{HP2)2{PH)2HPHP2H2P2 3.8307 



TABLE IX: Cfj^: 21-step ground state configurations common to all the energy-schemes, definitions are the 
same as the ones in TABLEHIII 



19 



U Sequence f 

Ui PHPzHz{PH)2{HP)3 3.9903 

U2 PHPzHz{PH)2{HP)z 4.2740 

U-i PHP3Hz{PH)2{HP)z 4.4030 

Ui PHPsHa{PH)2{HP)a 4.0885 

TABLE X: Designing sequence for 17-step configurations, Cl^^^: 12132121341312313 and Cl^'': 
12142121431412414. f: maximum gap-parameter. Subscript of H or P denotes the number of times the corre- 
sponding symbol is repeated. 



u 


E.AU) 


f 


Ui 


-11 


5.1636 


U2 


-80 


4.8322 


f/3 


-84 


5.1958 


U2 


-76 


4.3359 



TABLE XI: Square lattice walks. Number of contacts = 16. Common designing sequence is 
H2P4H2PHP2H2PHP2H2PHP2H; designable ground state configuration: 112234323341434412141123. 
Eg{U): minimum energy for the scheme U ; f: maximum value of the gap-parameter. 
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FIG. 6: (a) The ground state configuration 112234323341434412141123 displaying the designed sequence, 
H2P4H2PHP2H2PHP2H2PHP2H {H: red; P: blue), (b) Specific Heat as a function of temperature. 
Insets show the zeros of the partition function Ea.(ll4t for the energy schemes U2, U3 and U4. The sharp peak 
at T ~ 0.4 (black) corresponds to C/i . A single broad peak at T 2.5 is seen for both U2 (red) and C/4 (green). 
In contrast, data for Us (blue) show a prominent peak at T ~ 1.5 besides a broad peak at T ~ 3 (indicated by 
arrows). 
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FIG. 7: The designed sequence is HP2HPH2P2HPH2P2HPH2P2HPH2 (H: red; P: blue), 
(a) The ground state configuration 112233321111444433332211 displaying the designed sequence, is 
only for the scheme f/4 with a designability of 195 sequences. (b) The ground state configuration 
1 12234323341434412141 123 for all the energy schemes. While this sequence is the only sequence identifying 
this configuration as the ground state for schemes f/i — (73, it turned out to have the maximum gap-parameter 
(= 4.4151) for the scheme Ua- 



